Effect of weak disorder on the ground state of uniaxial dipolar 
spin systems in the upper critical dimension 
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Extensive Monte Carlo simulations are used to investigate the stability of the ferromagnetic ground 
state in three-dimensional systems of Ising dipoles with added quenched disorder. These systems 
model the collective ferromagnetic order observed in various systems with dipolar long-range in- 
teractions. The uniaxial dipolar spins are arranged on a face-centred cubic lattice with periodic 
boundary conditions. Finite-size scaling relations for the pure dipolar ferromagnetic system are 
derived by a renormalisation group calculation. These functions include logarithmic corrections to 
the expected mean field behaviour since the system is in its upper critical dimension. Scaled data 
confirm the validity of the finite-size scaling description and results are compared with subsequent 
analysis of weakly disordered systems. A disorder-temperature phase diagram displays the preser- 
vation of the ferromagnetic ground state with the addition of small amounts of disorder, suggesting 
the irrelevance of weak disorder in these systems. 
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INTRODUCTION 

Certain closely-packed configurations of uniaxial 
dipoles can undergo a transition to a ferromagnetic 
ground state even in the absence of a direct exchange in- 
teraction |l| . Dense arrangements of ferromagnetic mon- 
odomain particles and some molecular magnetic crystals 
are exemplary cases of this phenomenon [2(. However, 
the existence of this dipolar ferromagnetic (DFM) ground 
state is highly dependent on the arrangement of spins, 
due to the directional dependence of the dipole inter- 
action. Dipoles distributed on face-centred and body- 
centred cubic lattices [3|, |4| are known to exhibit the tran- 
sition to ferromagnetic ordering and it is straightforward 
(albeit computationally expensive) to simulate the tran- 
sition using the Monte Carlo method. Simulation results 
show a tendency to order at a well-defined temperature. 
Such simulation models also exhibit rapid deterioration 
of the DFM ground state with the introduction of disor- 
der, as shown in Fig.^ In this example, quenched disor- 
der is added to a closely-packed system of Ising dipoles 
in the form of randomly-signed exchange between near- 
est neighbours, as in the short-range Edwards- Anderson 
(EA) spin glass 5]. Aside from the limit in which these 
additional interactions dominate, the behaviour leading 
to the destruction of the ferromagnetic order is unknown. 

Harris 6] proposed a general criterion for changes to 
the behaviour of pure systems with the introduction of 
randomness. A new type of critical behaviour is pre- 
dicted, if the transition observed in the pure case is sharp 
in a specific sense. The transition must display a power- 
law divergence of the specific heat (with exponent a > 0), 
as opposed to a singular thermodynamic behaviour (with 
a < 0), where the transition is signalled only by a cusp 
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FIG. 1: Temperature dependence of magnetisation measure- 
ments for a FCC dipolar system with N = 108 spins for vary- 
ing amounts of disorder, controlled by the coupling constant, 
Jea, which describes the strength of the randomly-signed 
short-range exchange. 



in specific heat. Chayes et al. [7j examine a host of sys- 
tems for which a positive specific heat exponent in the 
pure system corresponds to a negative exponent in the 
random equivalent. Indeed, if the converse is true, and 
the pure system has a negative specific heat exponent, 
then the addition of disorder simply leads to irrelevant 
corrections to the scaling behaviour ||. 

Three-dimensional Ising dipolar systems have been 
shown to behave in a way characteristic of a mean field 
description 0, 0, 0] , whereas those comprising vector 
dipolar spins typically resemble systems interacting via 
short-range exchange [l^. This discrepancy stems from 
the fact that the upper critical dimension (defining the 
lower limit for which fluctuations can be neglected) is 
three in the former case and four in the latter. The 
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Ising dipolar system therefore represents a unique case 
for which predictions based on exact solutions from mean 
field theory at marginal dimensionality can be tested 
experimentally Furthermore, it presents an ideal 

testing ground for the extent to which the interaction 
range influences critical behaviour 0|. The marginality 
of the system manifests itself in the appearance of loga- 
rithmic corrections to mean field theory which produce 
a logarithmic divergence in the specific heat, implying 
a sharp transition, consistent with experimental findings 
for three-dimensional systems of Ising dipoles |l5| . 

This investigation examines the effect of disorder on 
a three-dimensional system of Ising dipoles by means of 
Monte Carlo simulations. Finite-size scaling is employed 
to locate the temperature of the expected phase transi- 
tion in the pure system and characterise its critical prop- 
erties. Data from weakly disordered systems are scaled 
in a similar manner in order to gauge the relevance of the 
disorder. It is clear, both from the raw data and from the 
scaling analysis presented, that the ferromagnetic ground 
state cannot be detected beyond a certain disorder level. 



MONTE CARLO SIMULATIONS 



The model Hamiltonian for the uniaxial dipole system 
with the Ising spin axis oriented in the z direction 
is, 
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(1) 

The vector S has length a and runs over c nearest neigh- 
bours, coupling spins via a quenched random bond vari- 
able, a r s £ {— 1,+1}. The magnitude of the exchange 
energy is denoted Jea and G = \{g^B) is the strength 
of the dipole interaction. Spins are placed on a lattice 
array with face-centred cubic (FCC) structure and num- 
ber density of one spin per unit cube, in the length scale 
of the simulated system. The dipole axes are aligned 
with one of the edges of the cubic lattice cells. Periodic 
boundary conditions are implemented by Ewald summa- 
tion for a cubic simulation box [Tt| . Energy, temperature 
and Jea are all measured in units of dipolar interaction 
energy for a pair of spins of unit length separation. 

Using the Metropolis Monte Carlo method, data are 
obtained for systems of linear dimension L = 3, 4, 5, 
6 and 7 via a process of slow-cooling. This is executed 
through a series of isothermal runs, in temperature steps 
of 0.025 and 0.05. After an initial prologue for equilibra- 
tion, several thermodynamic properties are measured, in- 
cluding magnetisation, (to); second and fourth magnetic 
moments, (to 2 ) and (to 4 ); energy and its square, (E) 
and (E 2 ); susceptibility, \ — (( m2 ) — ( TO ) )/T\ specific 
heat, c = ((E 2 ) — (E) 2 )/T; and the Binder cumulant, 



[7 = 1 — (to 4 ) /3(to 2 ) . Disorder is introduced in the 
form of quenched random interactions between nearest 
neighbours. These are short-range E A- type spin glass 
interactions of variable strength, Jea, and disorder av- 
erages are determined from tens to hundreds of disorder 
realisations. 



RENORMALISATION GROUP 
TRANSFORMATION 

The critical properties of the system are determined 
by scaling the data to achieve size-independent collapse 
of the temperature dependence of all observables. The 
required finite-size scaling Ansatze are derived by exam- 
ining the behaviour of the system under renormalisation 
group transformation. Since the test system comprises 
purely dipolar interactions, the random variable is omit- 
ted from the Hamiltonian in Eq. QJ. Whilst the frame- 
work for this transformation was developed many years 
ago 0| , the form of the scaling functions for the marginal 
case of the Ising dipolar system is yet to be explicitly 
derived. One can formalise the 'phenomenological renor- 
malisation' 0] occuring at the transition by transform- 
ing the lattice size L as L' = L/b. 

The modified Hamiltonian is converted to one describ- 
ing continuous spins via the standard <p A theory in which 
a coefficient uq acts as a control parameter for the in- 
teraction [l^. A source term, h, is included to enable 
functional differentiation of the free energy. The result 
takes the form, 
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where the following integral convention is employed here- 
after, J — > J <o<A i?J\ d w ith A ^ 7r/a. The function in 



iq Jo<g<A 

the first integral of Eq. @ is defined as 
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where the coefficient ro oc t = (T — T c ) /T c and T c is the 
critical temperature of the DFM transition. The remain- 
ing coefficients arise from an expansion of the dipolar 
term in the Hamiltonian in the limit of small q. The spin- 
spin correlation function for the Gaussian approximation 
to this model is simply the inverse of this function, Uo(q). 
For q — 0, the correlation function is shape-dependent 
and takes the form, 



G o (0) - (r + gold- 1 - (d - 2)*/oi))" 



where, 
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and 6ij is the angle between and the z axis 10]. The 
renormalisation is undertaken by integrating out large 
wavevector modes and expanding the non-Gaussian part 
of the Hamiltonian perturbatively for small Uq [2jj ■ The 
result is rescaled according to q = q' /b and <fr — C,<j)' qll 
where £ is chosen in order to keep the coefficient of the 
q 2 term constant to leading order in u 2 . The renormali- 
sation group equations can be written in differential form 
by taking b = e l and building up infinitesimal transfor- 
mations using the group composition property. 

The resulting scaling function for the free energy den- 
sity takes the form, 

/(*, h, u, L- 1 ) = b- d f(t', h', u', b L- 1 ) + g(t, h), (6) 

where g represents the analytic part of the transforma- 
tion and coefficients r$ and uo have been identified with 
variables t and u. The latter is a dangerous irrelevant 
variable since the free energy can become singular below 
the transition if it is removed. The correct finite-size scal- 
ing properties in this case are obtained by renormalising 
the system size to unity [ilj. A further substitution of 
4>' = (fi/u' 1 ^ 4 then gives a new scaling form for the free 
energy, 



f(t',h',u',l) + g(t,h)=f(t,h), 



i=t'/u' 1/2 ,h 
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The analytic part of the transformation also contributes 
to the singular dependence of the free energy on t and 
this is absorbed into the function /. In order to write 
the scaling forms in a representation conducive to data 
fitting, it is convenient to absorb all constant factors into 
the singular function and define a fitting parameter, v, 
which controls the shift imposed on the reduced temper- 
ature variable, t j2^. Combining Eq. JfJJl and Eq. J7J|, 
the free energy includes logarithmic corrections and is 
expressed as, 



/(t, h, u, L- 1 ) = L~ d f(t, h) + g(t, h), 



(9) 



t = L yt log 1/6 L(t + vL~ Vt log~ 2/3 L),h — L Vh h\og 1/4 L, 

(10) 

where yt = 3/2 and yh = 9/4 define the renormalisation 
exponents for the mean field theory and logarithmic func- 
tions are written in terms of some unitary length scale 

El. 



FINITE-SIZE SCALING 

The scaling forms for physical quantities follow then 
from simple derivatives of this function. They can be 



expanded in Taylor series about the finite-size limit in 
the critical region, given by small t and finite L [23| . The 
functions are, 
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for the Binder cumulant. The last term in the scaling 
form for the Binder cumulant arises due to the field de- 
pendence of the analytic part of the free energy. This cor- 
rection is not expected to affect the scaling of the other 
thermodynamic properties, as the logarithmic corrections 
will dominate their form [14|. Data collapse is achieved 
by plotting scaled observables with respect to the shifted 
temperature variable, t, so that data from systems on all 
length scales fall onto a single curve. 

The scaling form is only valid in the vicinity of the 
transition so the fitted data are restricted to this region 
in order to achieve meaningful collapse. The imposed 
restrictions are determined by pinning the set of each 
sample size to that of the smallest size. The subset of 
the smallest sample is chosen by minimising the weighted 
cost associated with fitting the restricted Binder cumu- 
lant data to the mean field approximation of Eq. (|14|) . 

The ability of a given set of fit parameters to ob- 
tain acceptable data collapse can be determined via a 
suitably robust ad hoc method. One scales every set 
of data points for all lattice sizes using the appropri- 
ate L-dependent logarithmic factors. Combining them 
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FIG. 2: (a) Finite-size dependence and (b) collapse of the Binder cumulant data in the critical region for the pure dipolar 
system. 



TABLE I: Parameter values for fitting the pure dipolar system, quoted with 99% confidence levels. 





T c 
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6i 
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Cl 


<m> 


3.15292 ±0.00013 


-1.5317 ±0.0011 


0.44 ±0.10 


-0.36 ±0.08 




<- 2 > 


3.15292 ±0.00005 


-1.14183 ± 0.00042 


0.79 ±0.17 


-0.62 ±0.13 




X 


3.152916 ±0.000004 


-1.14531 ± 0.00004 


0.08 ±0.02 


-0.07 ±0.01 




U 


3.152916 ±0.000005 


-1.64751 ±0.00004 


1.30 ±0.05 


-0.31 ±0.01 


-2.33 ±0.19 



in a single set, the data are then interpolated as a high- 
order polynomial defining an effective mean. The In- 
dependence of this function is implicit in the functional 
form of i =t(T,L;T c ,v). 

A cost function is formed from the error-weighted sum 
of squares for the vertical deviation of each scaled data 
point from this master fit curve. Minimising this func- 
tion numerically yields a set of parameters for each ob- 
servable. A global value for the critical temperature, T c , 
is obtained by averaging the weighted T c values from the 
initial fits and refitting all observables using the fixed 
average value. 



RESULTS 

In the case of the pure dipolar system, the dimension- 
less ratio U displays scale invariance at the transition 
with some statistical scatter. It is clear in Fig. Ela) that 
there is a size-dependent systematic shift away from a 
single scale invariant crossing at the critical tempera- 
ture. Fig. |2Jb) shows the collapse of these data, con- 
firming the validity of the derived scaling functions. The 
critical temperature and the fit parameter v are only ef- 
fective in collapsing the size-dependent data sets later- 
ally. The vertical displacement required to achieve to- 
tal collapse is controlled by parameters b\ and 62 of the 
logarithmic corrections and the importance of their in- 
clusion is vast when compared with previous analyses of 
Ising dipole systems |24|. Global fits for magnetisation 
and susceptibility also show excellent data collapse [2fj| . 



The correction arising from the analytic part of the free 
energy density transformation has little effect for all ob- 
servables aside from the Binder cumulant, and this is a 
further indication of the importance of the logarithmic 
corrections in the upper critical dimension of the uniax- 
ial dipolar system. It should be noted that the least sum 
of squares values corresponding to the parameters given 
by the global T c deviate by no more than 2% from their 
original values in the initial fits. 

Fitted parameter values are quoted with error esti- 
mates corresponding to 99% confidence in Table [I] The 
fitted values for the lateral parameters, T c and v, are 
calculated with much greater confidence than are those 
associated with vertical shift, which pertain to logarith- 
mic corrections. This is to be expected, since the refitting 
process is independent of T c , and the value chosen for v 
is pinned to that of the critical temperature. 

The mean field scaling functions yield good data col- 
lapse for a certain range of weak disorder. Fits are com- 
parable with those of the pure case up to a disorder 
strength of Jea ~ 0.25. As shown in Fig.|3{a), the mean 
field predictions are applicable at this value and the data 
sets collapse as they did in the pure case. Beyond this dis- 
order strength, the simulations continue to show a transi- 
tion to ferromagnetic ordering (see Fig.|3fb)). However, 
fitting data to the mean field scaling forms proves increas- 
ingly difficult in this intermediate range and confidence 
levels suggest the description is unsuitable here. When 
disorder exceeds a level corresponding to Jea ~ 0.5, 
there is no ferromagnetic transition and data cannot be 
scaled at all using the mean field predictions. This range 
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FIG. 3: (a) Collapse of the Binder cumulant data for the dipolar system with Jea = 0.25. (b) Temperature variation in Binder 
cumulant data for Jea ~ 0.40 suggests a ferromagnetic transition still exists. 



coincides with the disorder level at which the ferromag- 
netic transition appears to deteriorate in Fig. ^ 



DISCUSSION 

The critical temperature of the apparent ferromagnetic 
transition calculated by the fitting process decreases as 
the disorder strength is increased, as shown in Fig.0J The 
smooth change away from the value used to fit data for 
the pure case supports the preservation of the transition 
to DFM behaviour predicted by the mean field theory. 
Taking into account the numerical inaccuracy of the dis- 
ordered system data, this smoothness is evident in the 
case where added disorder does not exceed Jea — 0.25, 
whereupon there is an abrupt change from the low disor- 
der trend. This does not necessarily mean that systems 
beyond this disorder range are not of DFM character, 
since the mean field scaling functions can still achieve 
data collapse. However, discontinuities in the variation 
of parameter values may point towards a change in the 
critical behaviour. 

There are two plausible scenarios in which this change 
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could occur. The first predicts an intermediate transi- 
tion from pure DFM behaviour to disordered DFM be- 
haviour at some finite strength Jea — J*- One might 
suppose this new disordered phase exists in systems with 
disorder strength below some upper limit, Jea = J , at 
which spin glass ordering sets in. An alternative picture 
describes a gradual crossover to a new kind of disordered 
critical behaviour for any finite value of Jea- However, 
this scenario contradicts the results of scaling in the low 
disorder region to some extent. In particular, the onset 
of new behaviour with the introduction of disorder is not 
expected to be well-described by the mean field theory, as 
seen here. However, it is conceivable that the simulated 
systems are too small to detect such a subtle change. 

Whilst the behaviour of the system with intermediate 
disorder strength remains unclear, one can be confident 
that the strongly disordered systems simulated here do 
not have a DFM ground state. The behaviour observed in 
this region of phase space cannot be determined within 
the mean field description presented. Further analysis 
must address a revision of the scaling functions for the 
thermodynamic properties of disordered uniaxial dipolar 
systems, and an investigation of alternate forms of order- 
ing for the limit of strong disorder. 

We acknowledge the Australian Research Council and 
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